#include "sor.h"

void setPBoundaries(int imax,int jmax,double **P,int ** Flag, double dp);
void sor(
	double omg,
	double dx,
	double dy,
	int    imax,
	int    jmax,
	double **P,
	double **RS,
	double *res,
	int **Flag,
	double dp
)
{
	int i , j;
	double sor_calc , sum = 0;
	int count_Fluid = 0;
	
	sor_calc = 2.0 * ( (1.0 / pow(dx , 2) ) + ( 1.0 / pow(dy , 2 )));
	
	
	
	for(i = 1 ; i <= imax ; i++) {
		for(j = 1 ; j <= jmax; j++) {
			if ( is_fluid(Flag[i][j]) ) {
				
				P[i][j] = (1.0 - omg)*P[i][j] + (omg/sor_calc)*(((P[i+1][j]+P[i-1][j])/pow(dx,2))+((P[i][j+1]+P[i][j-1])/pow(dy,2)) - RS[i][j]);
				count_Fluid ++;
			}
		}
	}

	for(i = 1 ; i <= imax ; i++) {
		for(j = 1 ; j <= jmax; j++) {
			if ( is_fluid(Flag[i][j]) ) {

				sum += pow(((P[i-1][j] - 2*P[i][j] + P[i+1][j])/pow(dx,2))+((P[i-1][j] - 2*P[i][j]+P[i + 1][j])/pow(dx,2)) - RS[i][j],2);
			}
		}
	}

	*res = sqrt(sum/(count_Fluid));
}


void setPBoundaries(int imax,int jmax,double **P,int ** Flag, double dp)
{
	int i,j;
	/* boundary strip initialization */
		for(j = 1; j <= jmax ; j++) {
			if ( (Flag[0][j] == P_L) && (Flag[imax + 1][j] == P_R) ) {
				P[0][j] = 2*dp - P[1][j] ;
				P[imax + 1][j] = - P[imax][j];
			}
			else {
				P[0][j] = P[1][j];
				P[imax + 1][j] = P[imax][j];
			}
		}
		
		
		for(i = 1 ; i <= imax; i++) {
			P[i][0] = P[i][1];
			P[i][jmax + 1] = P[i][jmax];	
		}
		
		/* interior boundary cells initialization */
		for ( i = 1; i <= imax; i++ ) {
			for ( j = 1; j <= jmax; j++ ) {
				switch (Flag[i][j]) {
					case B_N:
						P[i][j] = P[i][j + 1];
						break;
					case B_W:
						P[i][j] = P[i - 1][j];
						break;
					case B_S:
						P[i][j] = P[i][j - 1];
						break;
					case B_E:
						P[i][j] = P[i + 1][j];
						break;
					case B_SW:
						P[i][j] = ( P[i - 1][j] + P[i][j - 1] ) / 2.0;
						break;
					case B_SE:
						P[i][j] = ( P[i + 1][j] + P[i][j - 1] ) / 2.0;
						break;
					case B_NE:
						P[i][j] = ( P[i][j + 1] + P[i + 1][j] ) / 2.0;
						break;
					case B_NW:
						P[i][j] = ( P[i - 1][j] + P[i][j + 1] ) / 2.0;
						break;
					default:
						break;
				}
			}
		}
	
}

